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Abstract 

We investigate the heat conduction properties of molecular junctions comprising anharmonic 
interactions. We find that nonlinear interactions can lead to novel phenomena: negative differential 
thermal conductance and heat rectification. Based on analytically solvable models we derive an 
expression for the heat current that clearly reflects the interplay between anharmonic interactions, 
strengths of coupling to the thermal reservoirs, and junction asymmetry. This expression indicates 
that negative differential thermal conductance shows up when the molecule is strongly coupled 
to the thermal baths, even in the absence of internal molecular nonlinearities. In contrast, diode 
like behavior is expected for a highly anharmonic molecule with an inherent structural asymmetry. 
Anharmonic interactions are also necessary for manifesting Fourier type transport. We briefly 
present an extension of our model system that can lead to this behavior. 

PACS numbers: 63.20.Ry, 44.10.+i, 05.60.-k, 66.70.+f 
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I. INTRODUCTION 



UnderstandinEf and controlling heat flow in nanoscale structures is of interest both from 
the fu„da.eMa. aspect ^ and fo. deviee appHcations S fl Q. The influential .ole of 

quantum effects and geometrical constrictions in low dimensional systems often results in 
fundamentally interesting behavior jsl . Recent theoretical and experimental studies demon- 
strated that the thermal transport properties of nanowires can be very different from the 
corresponding bulk properties P,!^- In the low temperature ballistic regime the phonon ther- 
mal conductance of a one-dimensional (ID) quantum wire is quantized, with g = rc'^k^T/Sh 
as the universal quantum conductance unit Q], where ks and h are the Boltzmann and 
Planck constants, respectively, and T is the temperature. Reflections from the boundaries 
and disorder in the wire can be further treated by considering a Landauer type expression 
for the heat current {h = 1) 0, 0, 13 

J = y duuT{u) [nL{uj) - nR{uj)] . (1) 

This relationship describes energy transfer between two (left (L) , right (i?)) thermal reser- 
voirs maintained at equilibrium with the temperatures Tl and Tr, respectively, in terms 
of the temperature independent transmission coefficient T{uj) for phonons of frequency uj. 
Here riKiyj) = (e'^^'^ — l) ^; Pk = ^/ksTK, {K = L,R), is the Bose-Einstein distribution 
characterizing the reservoirs. This expression assumes the absence of inelastic scattering 
processes, and the two opposite phonon flows of different temperatures are out of equilib- 
rium with each other. This leads to an anomalous transport of heat, where (classically) the 
energy flux is proportional to the temperature difference, AT = Tl — Tr, rather than to the 
temperature gradient, VT, as asserted by the Fourier law of conductivity, 

J = -ICAVT. (2) 



In this equation A is the cross section area normal to the direction of heat propagation and 
/C is the coefficient of thermal conductivity. An outstanding problem in statistical physics 
is to find out the necessary and sufficient conditions for attaining this normal (Fourier) 
law of heat conductivity in low dimensional systems jlll . 2, Q| • Among the 

crucial requirements explored is that the molecular potential energy should constitute strong 
anharmonic interactions. 
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Heat conductance experiments on short molecules or highly ordered structures provide 
results consistent with the Landauer expression. A micron length individual carbon nan- 
otube conducts heat ballistically without showing signatures of phonon-phonon scattering 
for temperatures up to 300 K 0, Q|. Intramolecular vibrational energy flow in bridged 
azulene-anthracene compounds could be explained by assuming ballistic energy transport in 
the chain connecting both chromophores jl9|. In contrast, calculations of heat flow through 
proteins show substantial contribution of anharmonic interactions leading to an enhancement 
of the energy current in comparison to the (artificial) purely harmonic situation [2^ . 

Anharmonic (nonlinear) interactions are also a tool for controlling heat flow in molecular 
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junctions with potential technological applications, e.g. a thermal diode 
thermal transistor Q]. We have recently demonstrated that when nonlinear interactions 
govern heat conduction, the heat current is asymmetric for forward and reversed temperature 
biases, provided the junction has some structural asvmmetr v |23| . 

In this paper we generalize the model developed in Ref. [2j|, and present a comprehen- 
sive analysis of the heat conduction properties of molecular junctions taking into account 
anharmonic interactions in the system. We discuss the influence of the following effects on 
heat flow through the junction: (i) interparticle potential, specifically the degree of molec- 
ular anharmonicity, (ii) molecule-thermal reservoirs contact interactions, and (iii) junction 
asymmetry with respect to the L and R ends. We derive an exact analytic expression for 
the heat current that clearly reflects the role of each of these factors in determining phonon 
dynamics. More specifically, we analyze the necessary conditions for demonstrating nega- 
tive differential thermal conductance (NDTC) and diode like behavior. We also follow the 
transition from the elastic Landauer formula to the Fourier law of conduction as anharmonic 
interactions are turned on. 

The paper is organized as follows: Section II presents our model system. Section III 
begins with a fully harmonic model and shows that it satisfies the Landauer formula. We 
then proceed and show that an asymmetric anharmonic molecule linearly coupled to thermal 
reservoirs, can rectify heat. Section IV further presents strong coupling models that exhibit 
NDTC. Under additional conditions, the nonlinear models extended to a I sites system satisfy 
the Fourier law of conductivity as explained in Section V. Section VI provides concluding 
remarks. 
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II. MODEL 



The model system consists a molecular unit connecting two thermal reservoirs left (L) and 
right (R) of inverse temperatures Pl — and I3r — T^^ respectively. Henceforth we take 
the Boltzmann constant as ks — 1- The general Hamiltonian includes three contributions: 
the molecular part (M), the two reservoirs {B) and the system- bath interaction {MB) 

H — Hm + Hb + Hmb- (3) 

For simplicity we assume that heat transfer is dominated by a specific single mode. The 
molecular term in the Hamiltonian is therefore given by 

Hm = ^ En\n){n\; En = nuJo, (4) 

n=0 

where uq is the frequency of the molecular oscillator {h = 1). We shall consider two situa- 
tions: harmonic model, and a two-level system (TLS) that simulates a highly anharmonic 
vibrational mode. For a harmonic molecule N is taken up to infinity. Strong anharmonicity 
is enforced by limiting n to 0, 1. The molecular mode is coupled either linearly (weakly) 
or nonlinearly (strongly) to the L and R thermal baths represented by sets of independent 
harmonic oscillators 

Hb^Hl + Hr; Hk^J2''^''W^ K^L,R. (5) 

Oj, Qj are boson creation and annihilation operators associated with the phonon modes of 
the harmonic baths. The L and R thermal baths are not coupled directly, only through their 
interaction with the molecular mode. We use the following model for the molecule-reservoirs 
interaction 

Af-l 

Hmb = Yl (^1^ ~ + - 1|) V^, (6) 

n=l 

where B are bath operators. This model assumes that transitions between molecular levels 
occur due to the environment excitations. Note that in general this interaction does not need 
to be additive in the thermal baths, i.e. we may consider situations in which B ^ B^ + Br, 
see Section IV. 
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The probabilities P„ to occupy the n state of the molecular oscillator satisfy the master 
equation 

Pn = - [nkd + {n + Pn 

+ (n + l)kdPn+i + nkuPn-i, (7) 

where the occupations are normalized Yin -^n = 1- Within second order perturbation theory 
the rates are given by 

/oo 
dre"^''^{B\r)B{0)), 
-oo 

/oo 
rfre-*'^°^(B(r)St(o)), (8) 
-oo 

where the average is done over the baths thermal distributions, irrespective of the fact that 
it may involve two distributions of different temperatures. In Section IV we demonstrate 
that these rates also apply in the strong molecule-baths interaction limit. 

A useful concept in the following discussion is the notion of an effective molecular tem- 
perature. It can be defined through the relative population of neighboring molecular levels 

= ~1 /p ° /TD \ - 

\0g{Pn+l/Pn) 

At steady state this ratio does not depend on n, see Eq. (|7j). We show below that the 
molecular temperature is given in terms of the reservoirs temperatures weighted by the 
molecule-baths coupling strengths. 

Given the reservoirs temperatures and Tr, we can define two other related parameters: 
the temperature difference AT = Tl — Tr and the average temperature Tq = (T^ + Tr)/2. 
The temperature difference can be experimentally imposed in various ways. Here we con- 
sider two situations: We may fix the temperature at the left reservoir while varying the 
temperature at the right side, 

(A) Tl = T, 

Tr = Ts- at. (10) 

For the same temperature difference we can also build a symmetric situation where the 
temperatures of both reservoirs are equally shifted, 

(B) Tl = Ts + AT/2 

Tr = T,- at/2. (11) 



The main difference between these two situations is that the average temperature is decreas- 
ing steadily with AT in the first case, while it is constant (T^) in (B). We will show below 
that these boundary conditions determine the effective mo/ecu/ar temperature which implies 
on the conduction properties of the system. 

Next we present the model Hamiltonians in the weak and strong molecule-bath interaction 
limits for either purely harmonic or a TLS molecular mode, and discuss the implications on 
the junction thermal conductance. 



III. WEAK SYSTEM-BATH COUPLING 



We begin by analyzing the heat conduction properties of a molecule coupled linearly to 
two thermal reservoirs of different temperatures |23|. The Hamiltonian is given by Eqs. 
0-® with linear (harmonic) system-bath interactions 

N-l 

Hmb = 5^ (fi|n-l)(n| + 5t|n)(n-l|) v^; 



n=l 



B = Bl + Br, (12) 
where the bath operators Bk satisfy 



Xj = {2ujj)~^''^{a] + aj)] K = L, R. (13) 

In the present linear coupling model, no correlations persist between the thermal baths and 
the rate constants (jH)) are additive in the L and R baths 

kd = kL + kn 

= fc^e-^^^" + fcije~^«"«, (14) 

with 

kK = Tk{'^q){1 + riKiuJo))] K = L,R. (15) 

Here n^(a;) = (e^^^'^-l)"\ Tk{uj) = Zljei^ ~ ^i) = ajV^mu^ Q, 

where m and c^o are the molecular oscillator mass and frequency respectively. 

The heat conduction properties of this model are obtained from the steady state solution 
of Eq. ((Tj) with the rates specified by Eqs. (fT H) -(fT3 |) . The steady-state heat flux calculated, 



e.g. at the right contact, is given by the sum 



N-l 



J = u;oJ2^ i^RPn - knP^-ie-"^''') , (16) 

n=l 

where positive sign indicates current going from left to right. 
A. Harmonic molecule 

For the harmonic model (A^ —>■ oo), putting = in (|7j), and searching a solution of 
the form P„ oc y"' we get a quadratic equation for y whose physically acceptable solution is 

y = 7— — T ' (17) 

which leads to the normalized state population 

P. = y"(l-1/). (18) 
Using Eq. (fTH|l we obtain the heat current (fTT)|) 

r r 

J = ujo ^ ^ (riL - nn) . (19) 
In the classical limit, uiq/Tk ^ 1 {K = L,R), it reduces to 

J=^^^{T,-Tr). (20) 

i L + J- _R 

This is a special case (with ^(V) = T r(T ^ + T r)^^ 6 {u — luq) consistent with our resonance 
energy transfer assumption) [2^ of the Landauer expression, Eq. ((T)). It is also consistent 
with the standard expression for the heat current through a perfect harmonic chain 2^. 

We emphasize on three important features of this result: (i) The heat current depends 
(classically) on the temperature difference between the two reservoirs, leading to divergent 
heat conductivity. Note that there is no need to introduce here the concept of the molecular 
temperature Tm- (ii) The current is the same when exchanging by F^, i.e. rectification 
cannot take place, (iii) The system cannot show the NDTC behavior, i.e. it is impossible 
to observe a decrease of the current with increasing temperature difference. This is true 
considering both models for the temperature drop- A and B, (fTUIl - lfTTI) . irrespective of the 
system symmetry. We can verify it by studying the AT derivative of the current 

dJ driL driR dn^ driji 



^ ^J\rTi a\rp a\rp + o/ \rp\- (^l) 



dAT dAT dAT OAT d{-AT) 



Since the term 



dAT Tl {e^^^o - 1)2 OAT '-^^^ 
is always positive (or zero), and similarly the second right hand side term, NDTC cannot 
show up in the fully harmonic model, and the heat current increases monotonically with the 
temperature difference. 

B. Anharmonic molecule 

We proceed to the case of a highly anharmonic molecule coupled -possible asymmetrically- 
but linearly, to two thermal reservoirs of different temperatures. We simulate strong anhar- 
monicity by modeling the molecular mode by a two levels system (TLS). The Hamiltonian 
for this model and the resulting rates are the same as presented throughout Eqs. ©-(HHI), 
except that we take n=0,l only. Following Eqs. ((Zl)-® we obtain the steady state levels 
population 

^i = r^' ^° = r4V- (23) 

We substitute it into Eq. ()16|) with N=2 and find the heat current 23| 

°ri(l + 2nL) + rR(l + 2n^) ^ ^ 

Next we calculate the molecular temperature Tm in the weak coupling-TLS case by sub- 
stituting the population into Eq. Q using Eqs. (fT^ - (fT^ . In the classical limit this 
results in 

TlTl + TrTr 

Tm = — . (25) 

L + i- R 

We can now study the implications of the different models for the temperature bias, Eqs. 
(jl(J|) - (jlip . on the conductance: In Model A the molecular temperature decreases monotoni- 
cally with the temperature difference 

Tif = T, - AT-^^. (26) 

In Model B we find 

T^B) , AT {Tl - Tr) 

J-M — -'s"! 7. 1 \^'} 

^ J- L + i _R 

which implies that for a sjTumetric (E/, = Yr) system, the molecular temperature is constant, 
whereas in the asymmetric situation it can either increase or decrease with AT, depending 
on the sign of Fj^ — Fjj. 
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In terms of the molecular temperature (Eq. ()25|) ). going into the classical limit, the heat 
current reduces into the simple form 

J=(T,-T„)^^^. (28) 

J- L + J- _R 

This relationship differs from the harmonic expression, ()2()|1 . by its implicit dependence on 
the internal molecular temperature. As we show next, this opens up the door for heat 
rectification and can also lead to the applicability of the Fourier law of conduction. If we 
still try to fit this expression into the Landauer form (^, we find that we have to define an 
effective temperature dependent trsnasmissioia coefficient T{uj,Tl,Tr) oc 1/Tm5{ijj — ojq). 

In Fig. n we display the current (Eq. ()24|l ) for a representative set of parameters. It 
increases monotonically with AT and saturates at high temperature gaps. We can verify 
this trend analytically as 



— 1 + 2nR - TTT^ 1 + 2ni 
[TL{l + 2nL) + Tn{l + 2nn)f 



dJ 
OAT 



which indicates that NDTC can not take place. However, Eq. ()24|) implies that the system 
can rectify heat current, i.e. the current can be different (in absolute values) when exchanging 
the reservoirs temperatures. Following [2^, defining the asymmetry parameter x such that 
= r(l -x);Tn = r(l + x) with -1 < X < 1 we get 



AJ = Jin = n;TR = Te) + J{Tl = T,-Tr = n) 



(30) 



(1 + + hrY - x'^i^L - npiY ' 
Here {T^) relates to the cold (hot) bath. Eq. (j3U|) implies that for small AT = Tl — Tr, AJ 
grows like AT^, and that the current is larger (in absolute value) when the cold bath is 
coupled more strongly to the molecular system. We exemplify this behavior at the inset of 
Fig.m 

We found therefore that a system consisting of an anharmonic molecular mode coupled 
linearly (harmonically) and asymmetrically to two thermal reservoirs of different temper- 
atures can rectify heat, though it cannot manifest the NDTC effect. NDTC requires an- 
harmonic interactions with the thermal baths, which may result in an effective nonlinear 
temperature dependent molecule-bath coupling term, see section IV. Therefore, there is no 
direct correspondence between these two phenomena. 




FIG. 1: Conduction properties of a TLS system in the weak coupling limit. ll'o=150 meV (full), 100 
meV (dashed), 25 meV (dotted). Ts=400 K (Model A), Tk=1.2 meV. Inset: Rectifying behavior 
of this model, uJo=25 meV, x=0.75 and rL=400 K, Tr = Tl - AT (full); TR=mO K, Tl = Tr-AT 
(dashed). 

C. General expression for the heat current 

We can generalize the harmonic ()20j) and anharmonic ()28p results and revise the current 
in the weak coupling limit (W) as 

Jw — uJo- — — — Js,B, [oi-) 

J- L + t /? J M 

where 

{1/2, Anharmonic TLS case ("Spin") 
(32) 
Tm/^o- Harmonic case ("Boson") 

For an intermediate anharmonicity we expect this function to attain an intermediate value, 
1/2 < fs^B < Tm /ujq. Note that fs^B can be retrieved by going into the classical limit of 
fs,B = [exp(co'o/^M) ± 1] Here the "spin" case takes the plus sign, and the "boson" situa- 
tion acquires the minus. It can be therefore interpreted as an effective molecular occupation 
factor. 

We can now clearly trace the influence of the different factors on the heat conductance. 
The thermal current is given by multiplying three terms: (i) A symmetric prefactor that 
includes the influence of the system-baths coupling, (ii) the factor ujq/Tm which includes 
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internal molecular properties: frequency and effective temperature, and (iii) the molecular 
occupation factor fs^B that varies between 1/2 for the strictly anharmonic system and Tm/ujq 
in the harmonic case. As we show next, the energy current has the same structure when 
system-bath interactions are strong. 



IV. STRONG SYSTEM-BATH COUPLING 



We turn now to the situation where the molecular mode is strongly coupled to the thermal 
reservoirs. As before, we discuss two limits, the harmonic case, and the anharmonic TLS 
situation. In both limits the model Hamiltonian includes the following terms, as in Eqs. 

N-l 

H = Y,En\n){n\ 

n=0 

+ V^Vn-i,n\n - l)(n|e^(""-^"-i) + c.c. 

n=l 

+ 5Z (33) 

jeL,R 

where En = ncuo, 17^ = 17^ + and ^1^ = i J2jeK ^n,j (^o] — j (K = L, R). In Appendix 
A we demonstrate that this model Hamiltonian equivalently represents a displaced molecular 
mode coupled nonlinearly to two thermal reservoirs. The coefficients A„ j are the effective 
system-bath interaction parameters that depend on the level index and the reservoir mode, 
see Appendix A. The Hamiltonian ()33p is similar to that defined in Eqs. Pj)- (fT^ . except that 
the L and R system-baths couplings appear as multiplicative rather than additive factors 
in the interaction term, implying non-separable transport at the two contacts [2^. The 
dynamics is still readily handled. For small V (the " non-adiabatic limit") the Hamiltonian 
leads again to the rate equation ((Tj) with 

h = \V\^C{ujo); K = \V\^Ci-cuo), (34) 

where 

/oo 
rfte^^«*C'(t), (35) 
-00 



11 



and 



C{t) = ^e«[f^"(*)~^n-i(*)]e-«[^"(o)-^^n-i(o)]^ 



R 



This may be evaluated explicitly to produce 



C{t) = CL{t)CR{ty, CK{t) = exp(-0^(t)), 



(36) 



(37) 



with 



(38) 



Note that we have omitted the n dependence from the rates above. This is supported by 
(i) taking all the inter-levels couplings to be equal, i.e. |V^_i,„| = V, and (ii) assuming that 
{Xnj — Xn-i,jY is the same for all n , e.g. A„j oc n, see Appendix A. 

Explicit expressions may be obtained using the short time approximation (valid for 
Sjgi<'(An,j — A„_ij)^ ^ 1 and/or at high temperature) whereupon 0(t) is expanded in 
powers of t keeping terms up to order t^. This leads to 



2tt 



exp 



2{Dl + Dl) 



(39) 



where 



E 



K 

M 



31 



Dl = Y,(\n,,-\n-l,,fuj]{2nK{uJ,) + l). 

jeK 



(40) 



E^ can be considered as the reorganization energy associated with the structural distortions 
of reservoirs modes around the isolated molecular vibration. In the classical limit {ujq/Tk 
0),DI = 2TkE^ 

Following Ref. 23| we calculate the steady state heat current utilizing 

N-l 



/oo 
duJiu[CRiiu)CLiuJo-Uj)Pn 

-Cr{-uj)Cl{-ujq + Uj)Pn-l]n, 



(41) 



12 



where 

C{uJo) = I dujCL{uJo - uj)Cr{uj), 



Cj,(uj) = p~("-^M)V4T^^iS (42) 

Eq. ()4H1 views the process |n) — > |n — 1) in which the molecular mode looses energy ujq as 
a combination of processes in which the system gives energy u (or gains it if < ) to the 
right bath and energy ujq — uj to the left one, with probability nCi{ujQ — uj)Cr{uj). A similar 
analysis applies to the process \n — 1) ^ \n). 

A. Harmonic molecule 

The levels population of an harmonic molecule (A^ oo) are calculated from the steady 
state solution of Eq. ((7j), leading to P„ = y"(l — y), y = k^/kd, with the rates conveyed by 
Eqs. pHl - lj^Oll . The heat current (jlT| is computed by first making the summation over n 

f-nP = 



oo 



Next, performing the integrals over frequency yields 

2^\V\'Ei,E^,iTL-Tn) 



J 



where 

= LMEij+E^j)/(Ej;jTL+E^,TR) _ l] ~\ (45) 

Before we discuss the heat conduction properties of this model we examine the anharmonic 
system. 
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FIG. 2: Conduction properties of the TLS system in the strong couphng hmit Tr = Tl — AT 
(model A), £^^=300 meV, {K = L,R) V=l meV, uJo=W meV. 

B. Anharmonic molecule 

The anharmonic model is described by the Hamiltonian ()33p with n = 0, 1. The steady 
state current is therefore obtained by reducing Eq. pTj) to 



J = \V\ 



duJLj[CR{uj)CL{uJo -uj)Pi 

oo 

-CRi-Uj)CLi-UJo + Uj)Po]. 



(46) 



Here Pq = C{ljq)/ {C{ljq) + C{—ljq)) and Pi = 1 — Pq are established from the steady state 
solution of (|7j) with the rates given by (jMj) . By following the same steps as for the harmonic 
system, the heat current (jlUj) is obtained as j23| 

2^\V\'Ej:,E^,iTL-Tn) 



J 



X e 



- (a;o - (Ef, ) ) V4{i?if Tfl) 



x/s 



(47) 



with the occupation factor 



fs 



'^MEj;,+Ef,)/{Ej;,Ti^+E^jTa) ^ ^ 



(48) 
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FIG. 3: Rectification in the strong coupling limit for a TLS system. £'Af=300 meV, V=l meV, 
a;o=10 meV, rL=300 K, Tr = Tl - AT (full); rR=300 K, Tl = Tr- AT (dashed). 

C. General expression for the heat current 

Next the harmonic (j44p and anharmonic results (j47p are reduced into a common form. We 
begin by evaluating the internal molecular temperature 0- In the present strong coupling 
case, for both harmonic and anharmonic molecular modes, it is given by 



e--/^- ^ P^^JP^ = (49) 
Using Eq. (jH^ we obtain the explicit expression 

2{El, + E^,) {Ei, + E^,) ■ ^'"^ 

The effective temperature in the strong coupling limit is therefore given by the algebraic 
average of the L and R temperatures weighted by the coupling strengths, here conveyed by 
the reservoirs reorganization energies. 

In terms of this quantity, we write a general expression for the current in the strong [S) 
coupling limit as 



= WW — ^-{u^o-Ei,-E^,)yATM{Ej;,+E^,) 

' ' ' iTM{Ei, + E^,f 
^ FL T 
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where 



(52) 



The plus sign relates to the anharmonic "spin" case, the minus stands for the harmonic 
"boson" situation. 

We analyze next the conduction properties of this model. Diode like behavior is expected 
when E^j ^ E^^, since then the resulting molecular temperature Tm is not the same when 
exchanging Tl with Tr. Note that in the present strong (nonlinear) coupling limit the 
molecule does not need to be strictly anharmonic for demonstrating this behavior, in contrast 
to the weak coupling situation. 

NDTC can also take place in the system, depending on the system asymmetry and the 
specific model for the applied temperature gradient. When the temperature bias is applied 
symmetrically at the L and R sides (model B, Eq. fill!) ). NDTC occurs for an asymmetric 
E^ 7^ E^j system. In model A the molecular temperature depends on AT even for a 
symmetric junction, providing NDTC. 

Figure |21 depicts an example of NDTC behavior in the system. The left reservoir is held 
at a constant temperature, while the temperature of the R reservoir is decreasing. We find 
that up to AT = Tl — Tn 100 K the current increases with the temperature bias, while 
above it, i.e. for lower Tr, the current goes down, and even diminishes (dotted line). 

We can also investigate the effect of asymmetrical contacts. We define the asymmetry 
parameter x such as E^;^ = Em{1 — x), E^j = -^^/(l + x)) < x < 1. Figure El presents the 
heat current when x 7^ 0- (^) Fo^' small x the current is almost the same for both forward 
and reversed operation modes, (b) At intermediate x values we find that for Tl = 100 K, 
Tr=300 K there is a maximal heat flow (dashed), while for the reversed operation (Tr = 100 
K, Tl=300 K) heat current is blocked (full), (c) For a highly asymmetric system heat flows 
predominantly in one direction. 

We can further formulate a general expression for the current that holds in both strong and 
weak interaction regimes and for either harmonic or anharmonic systems. For convenience, 
we copy here the weak {W) linear coupling result ()31|) 




fs,B- 



(53) 
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FIG. 4: (Color online) A schematic representation of the 1=4 molecular modes system with three 
intermediate thermal baths. 

Comparing it to Eq. (jKTjl guides us to the compact expression 

J = c^AT. (54) 

Here C includes the contact contribution, which is different in the weak and strong coupling 
regimes. It may depend on the molecule-baths microscopic couplings, molecular vibrational 
frequency and the reservoirs temperatures. It is not influenced by the degree of molecular 
harmonicity which affects only the "Spin-Boson" factor fs,B- The temperature Tm provides 
the effective temperature of the molecular system that is irrelevant in the fully harmonic 
case. We can therefore clearly distinguish in this expression between the role of the system 
harmonicity and the effect of molecule-baths interactions. 



V. FOURIER LAW OF CONDUCTION 

The validity of Fourier's law of heat conduction Q in ID lattices is an open issue jl^l ]. 
This law is a macroscopic consequence of ordinary diffusion at the microscopic level. Here we 
extend our single mode model and present a realization of a molecular chain that can lead to 
normal Fourier conduction, provided that the molecule is highly anharmonic, independently 
of the molecule-baths coupling strengths. We emphasize that we bypass the difficult task of 
showing how normal diffusion emerges in the system and simply assume non-correlated 
hopping motion between molecular units. Our sole mission here is to construct from the 
single mode result a model that supplies normal conduction. 

Figure E] depicts the system: We envision an array of / + 2 local heat baths where heat 
transfers along the chain of bath, single mode, bath, single mode...: L Bi —>■ B2 ^ ... ^ 
Bj R. The intermediate thermal baths Bi might be realized by large molecular groups 
where full thermalization to a local temperature occurs. This implies that anharmonic 
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interactions govern the dynamics within the intermediate baths. We also assume that the 
temperature difference Tl — Tr falls linearly on the system, i.e. the temperature at the i 
bath is 



Ti = Ts + 6T[i- {I + l)/2], AT = 5T{I + 1), (55) 

where i = (0, 1,2.../,/ + 1). The L and R bath are indiced by and / + 1 respectively, 
and 6T is the temperature difference between neighboring reservoirs. The resulting effective 
temperature of a molecular mode located in between each two baths z — 1 and i is 

= T, + 6T[2t - 1 - (/ + l)]/2, (56) 

assuming equal coupling along the chain (Fj = F or El/[ = Em)- 

The conductance of the unit /?j-i ^ Bi is defined through the relation J = Hi{Ti — 
Tj_i). In the diffusional hopping regime, the total conductance kj- is established by inversing 
the sum of all units resistances, kt = (X]j V'^*)"^- the linear coupling- TLS model (|HT|) 
it becomes 

where / = / + 1 is the actual length of the system, given by the number of internal molecular 
modes. This equation implies the validity of the Fourier law of heat conduction, J oc 
{Tl — Tr)/1. In addition, the prefactor depends on the inverse temperature as expected 
in the high temperature limit (2^ where interactions among phonons are dominant. The 
same result holds when the temperature falls mainly on the contacts whereas the temperature 
distribution along the central baths is almost constant. Strong molecule-bath interactions 
contribute a temperature dependent prefactor, but do not modify the form assuming 
AT/Ts ^ 1. A possible realization of this system is a fullerene polymer 



VI. CONCLUSIONS 



Using a simple theoretical model we have investigated the effect of anharmonic interac- 
tions on heat flow through molecular junctions. Our general expressions for the heat current 
(IHll) and (jHT|) clearly manifest the interplay between the system anharmonicity, system-bath 
coupling and junction asymmetry. We have also extended our model into a chain of / local 
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FIG. 5: An overview of the parameter ranges providing negative differential thermal conductance 
(NDTC) and diode like behavior. The x axis, (weak and strong coupling) relates to the system- 
bath interaction model. In order to obtain rectification the system should be asymmetric with 
respect to the L and R ends. 

molecular sites, and indicated that anharmonicity is a crucial requirement for achieving nor- 
mal conductance. We have found that nonlinear interactions can lead to novel phenomena: 
negative differential thermal conductance and heat rectification. NDTC takes place when 
the molecular mode is strongly coupled to the thermal environment. In contrast, diode like 
behavior originates from the combination of substantial molecular anharmonicities with a 
structural asymmetry. Figure El presents an overview of the different regimes studied, and 
the nonlinear effects observed in each case. 

We would also like to draw an analogy between the nonlinear behavior discussed in this 
paper and some nonlinear effects discovered in molecular level electron carrying systems: 
The negative differential resistance observed in molecular films of Ceo could be explained 
due to a voltage dependent tunneling barrier j^]. Rectification of electron current was 
theoretically exhibited in one-dimensional asymmetric electronic conductors with screened 
electron-electron interactions [sO^ . 

Control of heat flow through molecules by employing nonlinear interactions might be 
useful for different applications. In molecular electronic local heating of nanoscale devices 

n 

might cause structural instabilities undermining the junction integrity j^. Engineering 
good thermal contacts and cooling of the the junction are necessary for a stable and 
reliable operation mode. Control of vibrational energy transfer in molecules affects chemical 
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processes, e.g. reaction pathways, bond breaking processes, and folding dynamics j33|. Fi 



nally, we propose building technological devices based on heat flow, in analogy with electron 
current devices js^ . 
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APPENDIX A: Microscopic model for the strong coupling Hamiltonian 

The strong coupling Hamiltonian Eq. (|33|) can be derived from the following microscopic 
model 

N-1 J2. / N-1 ^2 

H = ^7Wo\n){n\ + 5^ Y a;^- - ^7^'^^^^' 

n=0 jeL,R \ n=0 J 

N-1 



+ VnVn,n-i\n) {n - 1| + VnVn-i^n\n - 1)( 
p 



n=l 

jeL,R 



which describes a forced oscillator of frequency ujq strongly interacting with the L and R 
thermal baths. The nonlinear contributions are concealed in the second element of ()A1|) 
providing high order terms such , with the molecular coordinate. Here Xj and 

Pj are the displacement and momentum of the reservoirs harmonic modes with frequency 
ujj, an,j is the system-bath coupling parameter and Vn,n-i is the effective inter-level matrix 
element. We can expand the quadratic term in ()A1|) and obtain 



N-1 / \ 

H = ''^^n j Uq — XjOnj 1 \n){n\ 

n=0 \ j€L,R J 



N-1 

+ ^ \/nVn,n-i\n) {n - 1| + VnVn-i,n\n -l){n\ 

n=l 

N-1 2 2 

ft fv 

+ J2u;,a]a, + J2Y1 ^1")^- (A2) 

jeL,R n=0 jeL,R J 
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iere Xj = (aj+aj) / ^J%j'j and p-,- = i\fIoJj2[a}^—aj). Use of the small polaron transformation 
a3, ^ = UHU-\ leads to 

n=0 n=0 j&L,R 3 

+ Yl VnVn-i,n\n - l)(n|e^(^"-^"-^) + c.c. 

n=l 

+ ^ cuj-ajaj, (A3) 
jeL,R 



where 



and where 



The term 



U = U':-,'a^, ?7„ = exp(-^a„|n)(n|), (A4) 



A„,,- = {2u^)-'/^nan,y (A5) 



N-l 2 2 

n=0 J J 

exactly cancels the oc term in Eq. ()A3|) . and we finally recover the strong coupling 
Hamiltonian ()33j) . 
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